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ABSTRACT 

We present a new method for the evaluation of the age and age-spread among pre-main-sequence 
(PMS) stars in star-forming regions in the Magellanic Clouds, accounting simultaneously for photomet- 
ric errors, unresolved binarity, differential extinction, stellar variability, accretion and crowding. The 
application of the method is performed with the statistical construction of synthetic color-magnitude 
diagrams using isochrones from two families of PMS evolutionary models. We convert each isochrone 
into 2D probability distributions of artificial PMS stars in the CMD by applying the aforementioned 
biases that dislocate these stars from their original CMD positions. A maximum-likelihood technique 
is then applied to derive the probability for each observed star to have a certain age, as well as the best 
age for the entire cluster. We apply our method to the photometric catalog of ~ 2 000 PMS stars in 
the young association LH 95 in the Large Magellanic Cloud, based on the deepest HST/ACS imaging 
ever performed toward this galaxy, with a detection limit of V ~ 28, corresponding to M ~ 0.2 M Q . 
We assume the Initial Mass Function and reddening distribution for the system, as they are previously 
derived by us. Our treatment shows that the age determination is very sensitive to the considered 
grid of evolutionary models and the assumed binary fraction. The age of LH 95 is found to vary from 
2.8 Myr to 4.4 Myr, depending on these factors. We evaluate the accuracy of our age estimation and 
we find that the method is fairly accurate in the PMS regime, while the precision of the measurement 
of the age is lower at higher luminosities. Our analysis allows us to disentangle a real age-spread from 
the apparent CMD-broadening caused by the physical and observational biases. We find that LH 95 
hosts an age-spread well represented by a gaussian distribution with a FWHM of the order of 2.8 Myr 
to 4.2 Myr depending on the model and binary fraction. We detect a dependence of the average age of 
the system with stellar mass. This dependence does not appear to have any physical meaning, being 
rather due to imperfections of the PMS evolutionary models, which tend to predict lower ages for the 
intermediate masses, and higher ages for low-mass stars. 

Subject headings: stars: formation — stars: pre-main sequence — Hertzsprung-Russell and C-M dia- 
grams — open clusters and associations: individual (LH 95) — methods: statistical 



1. INTRODUCTION 

The Large and Small Magellanic Clouds (LMC, SMC), 
the dwarf metal-poor companion-galaxies to our own, are 
well-suited environments for the study of resolved star 
formation at low metallicity. Indeed, both galaxies are 
characterized by an exceptional sample of star-forming 
regions, where hydrogen has been already started to be 
ionized by the UV- winds of massiv e stars, the Hn regions 



(Hcnize 1956 Davies et al. 1976). A plethora of young 



stellar associ ations and clusters are embedded in thes e 
regions (e.g., Bica & Schmitt 1995; Bica et al. 1999), 
giving evidence of recent (or current) clustered star for- 
mation. For many years, studies of the bright stellar 
content of these systems has provided and still provide 
unprecedented information on massive stellar evolution 
at low metallicities (see, e.g., Walborn et al.|1995| |2010| ). 
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On the other hand, their faint stellar content, which is 
still in the pre-main-sequence (PMS) phase, became only 
recently accessible exclusively th rough imaging with the 
Hubble Space Te lescope (see e.g., Gouliermis et al.|2006 
Nota et al.|[2006| ) 



The identification of the low-mass PMS stars of star- 
forming regions in the Milky Way is quite demanding, 
requiring various techniques such as variability surveys, 
spectroscopy, det ection of emission li nes and excess due 
to accretion (e.g., Briceno et al.|2 007), in order to decon- 
taminate the observed stellar samples from the evolved 
populations of the Galactic disk. This is not a prob- 
lem for the observational study of low mass population 
of young clusters in the Magellanic Clouds, where mem- 
berships can be assigned simply from the position of the 
stars in the color-magnitude diagram (CMD) alone. This 
is guaranteed by the negligible distance spread of the 
LMC or SMC field, in comparison to the distance of the 
two galaxies. 

Naturally, this observational advantage has important 
implications in our general understanding of clustered 
star formation and in consequence of PMS evolution, the 
Initial Mass Function (IMF) and dynamical evolution of 
young star-forming clusters. Clearly, in order to detect 
low-mass stars in the MCs, imaging with high resolution 
and high sensitivity in optical wavelengths is required, 
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and the Advanced Camera for Surveys (ACS) onboard 
the HST has been so far an ideal instrument for this 
purpose. Indeed, deep imaging at a detection limit of 
m 555 ~ 28 with ACS/WFC o f the young LMC cluster 
LH 95 (ILucke & Hodge 19701 proved for the first time 



that young MCs stellar sys tems host exceptional large 
numbers of faint PMS stars ( Gouliermis et al.|2007 ), pro- 
viding the necessary stellar numbers to address the issue 
of the l ow-mass IMF in star-forming regions of another 
galaxy (Da Rio, Gouliermis, & Henning 2009 1 . 



PMS stars remain in the contraction phase for sev- 
eral times 10 Myr, and therefore they are chronometers 
of star formation and its duration. However, the esti- 
mate of stellar ages in PMS regime from photometry is 
not without difficulties. This is because the observed 
fluxes are affected by several physical and observational 
effects. While these stars are optically visible, they are 
still surrounded by gas and dust from their formation, 
typically in circumstellar disks. Irregular accretion, due 
to instabilities in these disks, produces variability in the 
brightness of these stars. Moreover, PMS sources ex- 
hibit periodic fluctuations in light produced by rotating 
starspots. The variable nature of PMS stars dislocates 
them from their theoretical positions in the CMD. This 
effect in addition to observable biases, such as crowding 
and unresolved binarity, results in a broadening of the 
sequence of PMS stars in the CMD. As a consequence 
the measurement of the age of the host system and any 
intrinsic age-spread, which may indicate continuous star 
formation, cannot be derived from simple comparisons of 
PMS isochrones with the observed CMD. This problem is 
particulary critical for our understanding of the star for- 



mation processes in every region. In particular, Naylor 



(2009) showed that neglecting unresolved binaries may 
lead to an underestimation of the age of a young cluster 
of a factor 1.5-2. Concerning the star formation dura- 
tion, the apparent luminosity spread in the CMDs, if not 
disentangled from a possible intrinsic age spread, leads 
to an overestimation of the star formation timescale. 

With spectroscopic and variability surveys with HST 
in star-forming regions of the MCs being practically im- 
possible, the only way to disentangle the age and the 
duration of star formation in such regions from the rich 
samples of PMS stars is with the application of a statis- 
tical analysis, from a comparison between the observed 
CMDs and simulated ones. 

In this study we present a new maximum-likelihood 
method, for the derivation of the ages of PMS stars and 
the identification of any age-spread among them. We ap- 
ply our method to the rich sample of PMS stars detected 
with ACS in LH 95 to determine the age of the stars from 
their positions in the CMD in a probabilistic manner, ac- 
counting for the physical and observational biases that 
affect their colors and magnitudes. Our data set is pre- 
sented in § [2] and the evolutionary models to be used for 
our synthetic photometry are described in § [3j In order 
to achieve the most accurate results we further refine our 
previous selection of PMS stars members of the system 
in § |4j In § [5] we construct modeled 2D density distribu- 
tions of PMS stars in synthetic CMDs by applying the 
considered observational and physical biases to the evo- 
lutionary models, and we apply our maximum-likelihood 
technique for the derivation of the age of the system. In 
§ [6] we use our method to investigate the existence of an 



age-spread in LH 95, and we discuss our findings in §[7j 
The summary of this study is given in § [8j 

2. THE DATA SET 

In this study we make use of the photometry of the 
young association LH 95, obtained from deep Hubble 
Space Telescope (HST) imaging with the Wide-Field 
Channel of the Advanced Camera for Surveys (ACS) 
in the filters F555W and F814W, roughly correspond- 
ing to standard V- and /-bands respectively (program 
GO-10566; PI: D. Gouliermis). The reduc tion of the im- 
ages and their photometry is discussed in IGouliermis et 
al. (|2007[), and is presented in its entirety by Da Rio 
Gouliermis, &: Henning (2009 from hereafter Paper 1). 
T he photometric completeness reaches the limit of 50% 
at 771555 — 27.75 mag and mgi4 ~ 27 mag, reduced due 
to crowding by ~ 0.5 mag in the central area occupied 
by the association. A control field ~ 2' away from the 
system was also observed with an identical configuration. 
It has been used in Paper I to isolate the true stellar pop- 
ulation of the association from that of the general LMC 
field by performing a statistical field subtraction on the 
observed color-magnitude diagram (CMD). The HST ob- 
servations of the two fields and in the two bands were not 
consecutive. In particular, the /—band imaging of LH95 
was carried out 4 days after the V— band. This will be 
crucial in the analysis we present in this work to model 
the observed effect of rotational variability. 

The derived cluster-member population includes 2 104 
stars located in the PMS and Upper main-sequence 
(UMS) region of the CMD with corresponding masses as 
low as ~ 0.2 Mq. In Paper I we isolate the central area 
of the observed field-of-view (FoV) , which covers a sur- 
face of 1.43 arcmin 2 (~ 260 pc 2 ) corresponding to about 
12.6% of the whole FoV, as the most representative of the 
association LH 95. This area, which encompass the vast 
majority of the observed PMS population, is referred to 
as the "central region' 7 . 

Given low density of stars in LH95 (of the order of 10 
stars/pc 2 down to 0.3 — 0.2 M ), it can be considered 
a loose association. In our previous work three denser 
substructures - referred to as "subcluster A, B, C" - were 
identified. They include about half of the stellar popula- 
tion, while another half is distributed within the central 
region. In Paper I we roughly estimated the dynami- 
cal properties of the central region as well as the three 
subclusters. This was done estimating the total stellar 
mass, extrapolating our IMF for masses below our detec- 
tion limit, and assuming a velocity distribution for the 
stars of the order of ~ 2 km/s. For the entire region, the 
crossing time is comparable with the age of the system, 
and it is few times shorter when computed singularly for 
each of the three subclusters. It was still unclear whether 
these stellar agglomerates were formed independently on 
each other, from separate star forming events, or if their 
spatial distribution was simply caused by a clumpy dis- 
tribution of the ISM in the parental cloud, during a single 
star formation event. In Section |5.3| we investigate this, 
excluding the presence of a spatial variability of the stel- 
lar ages. 

3. DESCRIPTION OF THE EVOLUTIONARY MODELS 

For the purposes of our study we calculate a set 
of evolutionary models for PMS stars with the use of 
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Fig. 1. — A sample of isochrones computed by us using the 
FRANEC code. The shown isochrones correspond to ages of 0.5, 
1, 2, 3, 5, 10, 20, 50, 100 and 200 Myr. 



the FRANEC (Frascati Rap hson Newton Evolutionary 



Code) stellar evolution cod e (Chieffi & Straniero 1989 
i Degrinnocenti e~ al. 2008). Opacity tables are t aker 
from Ferguson et al. (2005) for log T[K] < 4.5 and Igle 
|sias 
tion of 



Rogers (1996) for higher tem peratures. The equa - 
t state (EOS) is described in|Rogers et al. (1996) 



Both opacity tables and EOS are calculated to r a heavy 
elements mix ture equal to the solar mixture of |Asplund| 
et al. (2005). Our models are then completely self- 



consistent, with an unique chemical composition in the 
whole structure. 

The value of the mixing length parameter, a, has been 
chosen by fitting a solar model, obtained to be a = 1.9. 
Helium and metal content of the models are those best 
representative for the LMC, i.e. Y = 0.25 and Z = 0.01 
respectively. The evolutionary tracks have been calcu- 
lated for a range of masses between 0.2 and 6.5 M . We 
have interpolated a dense grid of isochrones from the evo- 
lutionary tracks, which covers ages from 0.5 to 200 Myr, 
with steps of 0.1 Myr below 10 Myr, with steps of 0.5 Myr 
between 10 and 20 Myr, with steps of 1 Myr between 20 
and 50 Myr isochrones, and with steps of 5 Myr above 
50 Myr. For ages r < 1.6 Myr the isochrones do not 
reach the main sequence (MS) for the highest masses in 
our sample, while at r ~ 50 Myr the models cover also 
post-MS evolutionary phases for intermediate masses. 

The evolution of light elements (Deuterium, Lithium, 
Beryllium and Boron) is explicitly included in this ver- 
sion of FRANEC code. Deuterium burning plays an 
important role in the very early stages of PMS evolution 
for stars with M < 0.6 M . For these stars the energy 
provided by the D(p, 7) 3 He reaction is able to slow down 
the contraction in the initial part of the Hayashi track. 
This can be seen in Fig. [II where it can be seen that the 
low-mass end of the PMS isochrones are bent towards 
higher luminosities a for ages younger than ~ 5 Myr. 

In order to explore the dependance of our findings 
on the assumed models, the consideration apart from 
FRANEC o f additional families of PMS evolutionary 
models (e.g., i D'Antona 
Stahlcr 



1997 Palla 



Mazzitelli 1994 



1 



1999 



Baraffe et al. 
would be justified. How- 
ever, the majority of the available grids of models are not 



available for low metallically populations like those in the 
LM C. Th e grid of PMS evolutionary models by |Siess et 
al. ( 2000 ) is the only one that has been computed also for 
a subsolar value of metallicity Z = 0.01 and covers the 
mass range which is interesting to us, from the hydro- 
gen burning limit of ~ 0.1 M to the intermediate-mass 
regime with M 7 M . As a consequence, in ou r sub- 
sequent ana lysis we consider both FRANEC and Siess 



et al. (|2000[) models. We will refer to the latter simply 



as the Siess models. 

3.1. Synthetic Photometry 

In order to compare our photometry with the evolu- 
tionary models, we convert them from the physical units 
(Teg, log Lboi) into absolute magnitudes in the HST/ACS 
photometric system, by means of synthetic photometry 
on a grid of atmosphere models. As discussed thoroughly 
in Paper I, the surface gravity dependence in the compu- 
tation of optical colors for PMS stars must be accounted 
for, especially in the low-mass regime. Since atmosphere 
models arc generally expressed in terms of T c ff and log 
synthetic photometry provides a precise consideration 
of the color dependence on surface gravity, allowing to 
isolate one unique spectrum as the best representative 
of the stellar parameters for each point in a theoretical 
isochrone. 

The choi ce of a reliable atmo sphere grid is critical for 
cool stars (|Da Rio et al. 2010), whose spectral energy 



distributions (SEDs) are dominated by broad molecular 
bands. Moreover, these spectral features are gravity de- 
pendant, and since the stellar surface gravity varies dur- 
ing PMS contraction, optical colors are age dependant for 
young stars. In Paper I for the d erivation of the IMF of 
LH 95 we utilized the NextGen flHauschildt et al.|1999 1 
synthetic spectra to convert the theoretical Siess models 
into c olors and magnitudes, extended with the |Kurucz| 
(1993) grid for the highest temperatures (T eff > 8000 K). 
While our ACS photometry of LH 95 in Paper I includes 
PMS stars down to ~ 0.2 M , the detection incomplete- 
ness could not be corrected for masses below ~ 0.4 M . 
Therefore, the IMF derived in Paper I was constructed 
for masses down to ~ 0.4 M . However, in the present 
study on the age of LH 95 and the age distribution of 
its low-mass PMS stars, the whole sample of detected 
stars down to ~ 0.2 M can be used. Therefore we want 
to refine the atmosphere models we adopt to compute 
i ntrinsic colors for the very low stellar masses. 



Baraffe et al. ( 1998 1 show that NextGen models can- 
not reproduce empirical solar-metallicity V vs (V — I) 
color-magnitude diagrams for main sequence popula- 
tions. The problems set in for (V — I) > 2, which corre- 
sponds to T c g < 3600 K, and can be attributed to short- 
comings in the trea tment of molecules. In particular, 



Allard et al. (2000) illustrate that the main reason for 
the failure of NextGen spectra in matching the opti- 
cal properties of late-type stars is the incompleteness of 
the opacities of TiO lines and water. In that work they 
present new families of synthetic spectra (AMES mod- 
els) with updated opacities, and discuss the comparison 
with optical and NIR data, finding a good match with 
e mpirical MS dw arf colors from the photometric catalog 
o f|Leggett|(|1992D. 



|Da Rio et al. 12010 ) , in their analysis of a multi-color 
optical survey of the Orion Nebula Cluster down to the 
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hydrogen burning limit, show that, for young, M-type 
stars: a) the optical colors are bluer than those of MS 
dwarfs at the same temperature; b) intrinsic colors ob- 
tained from synthetic photometry using AMES spec- 
tra better reproduce the observed colors; c) for late M- 
type stars the AMES colors are even slightly bluer than 
the observed sequence, needing an additional correction 
based on the empirical data. As we show later on, the 
predicted colors from NextGen models are even bluer 
than the AMES ones for young stars. Mo reover, the 
corr ection to the AMES colors proposed in |Da Rio et] 
al. (2010) is critical for temperatures below our detec- 
tion limit. As a consequence, we are confident that for 
the analysis in this paper, theAMES models are in any 
case more appropriate than the NextGen, and their still 
present uncertainties do not play a significant role when 
used on our LH 05 photometry. 

Unfortunately the AMES library is available only for 
solar metallicity, and therefore we proceed by perform- 
ing an approximation that provides the best estimate of 
the correct magnitudes and colors for M-type stars with 
metallicity Z = 0.01. In conclusion, we perform syn- 
thetic photometry on the evolutionary models considered 
here b y applying the technique presented by|Da Rio et al.| 



(2010) by making use of three families of synthetic spec- 
traT(a) NextGen, for a LMC metallicity of Z ~ 0.01, 
(b) NextGen for solar metallicity Z = 0.02, and (c) 
combined NextGen for T cff > 5000 K and AMES for 
T cff < 5000 K for solar metallicity. For T fiff > 8000 K 
all three used grids are extended using Kurucz ( 1993 1 
spectra. 

Results of our transformations are shown in Fig.^/e/f) 
for a number of FRANEC isochrones with ages spanning 
from 0.5 to 100 Myr. In this figure it is shown that while 
the dependence of colors on the used synthetic spectra 
is negligible for low color-terms that correspond to high 
T e ff, the offset increases for late spectral types and young 
ages. In particular, the more "correct" AMES atmo- 
spheres are up to 0.2 mag redder than the NextGen 
ones for the same metallicity, and a decrease of [M/H] 
by a factor of two (from the solar to LMC metallicity) 
leads to a similar offset towards bluer colors. 

In order to estimate the magnitudes and colors of 
AMES-transformed isochrones for the LMC metallic- 
ity, we compute the displacement produced by chang- 
ing the metallicity from Z = 0.02 to Z = 0.01 with the 
use of the NextGen spectra for every point of every 
isochrone. These displacements are represented with ar- 
rows in Fig. J2jrig/i£/ They can be considered as a fairly 
good approximation of the real change in the observed 
colors and magnitudes in relation to T c g and log g, caused 
by the decrease of metallicity. Therefore, we apply the 
same offset, point by point, to isochrones transformed 
with the use of AMES spectra with solar metallicity, in 
order to derive their counterpart isochrones for Z = 0.01. 

Fig. [3] shows a sample of Franec isochrones trans- 
formed to the observable plane with the use of the 
AMES spectra library extrapolated to the LMC metal- 
licity plotted over the CMD of LH 95. A distance mod- 
ulus of m — M — 18.41 mag and average extinction of 
Ay = 0.6 mag, derived in Paper I, are applied to the 
isochrones. The PMS stars of LH 95, roughly following 
the 5 Myr isochrone, exhibit a wide apparent spread that 
spans from the youngest ages available in the models to 



the main sequence. 

4. FIELD SUBTRACTION REFINEMENT 

Our observations of the control field 2' away from 
LH 95 with setup identical to those of the system al- 
lowed us in Paper I to perform a thorough statistical 
removal of the contamination of the observed population 
by the stars of the local LMC field, on the verified basis 
that there are no PMS stars in the LMC field. Indeed, 
this method provided a sufficiently clean CMD of the 
PMS cluster-members of LH 95 (shown in Fig. [3]) . How- 
ever, there is an anomalous presence of a smalTnumber 
of sources in the field-subtracted CMD located between 
the general PMS population and the MS. These sources 
can be seen in Fig. [3] located in the neighborhood of the 
100 and 30 Myr isochrones. In Paper I we hypothesized 
the origin of these contaminants as a binary sequence of 
the field MS population, and we demonstrated that their 
removal does not affect the shape of the derived IMF. 

However, for the present study, which focuses on the 
identification of any age distribution among PMS stars, 
the nature of these sources and their membership should 
be reevaluated. In order to study the spread of the field 
MS in the CMD of both LH 95 and the control field, 
we construct the Hess diagram of the complete stellar 
sample observed in the area of the association before 
any field subtraction. We bin the stars in the color- 
and magnitude-axis of the CMD and smooth the counts 
with a Gaussian kernel of 0.05 mag in 771555 — msii and 
0.16 mag in m 814 . We repeat this process to construct 
the Hess diagram of the control field. The constructed 
Hess diagrams are shown in Fig. [4j An important dif- 
ference between the field population covered in the area 
of LH 95 and that in the control field is that LH 95 is 
embedded in th e H11 region n amed LHA 120- N 64C (or 
in short N 64C, |Henize1[1956| or DEM L 252 ( |Davies et 
al.|1976 ) , and therefore the reddening of the background 
field stars in the area of the system should be somewhat 
higher than those in the control field. However, a com- 
parison of the Hess diagrams of Fig. [4] shows that there 
is no significant variation of the field population between 
the two areas, and therefore the appearance of the neb- 
ula does not seem to affect significantly the field MS. 
This suggests that the LMC field population is mainly 
distributed in foreground with respect to LH95. Nev- 
ertheless, a small fraction of field stars could be located 
behind the nebula, and be affected by a larger extinction. 

In order to test this possibility and further eliminate 
the sample of PMS cluster-members from potential field 
MS contaminants, we consider diagonal cross-sections 
of width 0.33 mag in the direction defined by the line 
m8i4 = —1.92 x (777555 — TO814), and we fit with Gaussians 
the distributions of the stars along each of these strips. 
A sum of two Gaussian distributions is considered for 
the fit of the CMD of LH 95 and a single Gaussian for 
that of the control field. The centers and standard de- 
viations, cr, of the best-fitting Gaussians are drawn over 
the Hess diagrams of Fig. H] We compare the a of the 
best-fitting Gaussians to the MS, as a measurement of 
the broadening of the MS, between the CMD of LH 95 
and that of the control field in Fig. [5] where we plot a 
for each selected strip as a function ofthe corresponding 
average luminosity. In this figure can be seen that a and 
therefore the MS broadening is systematically smaller in 
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Fig. 2. — Left panel: CMD showing our FRANEC isochrones for ages of 0.5, 1, 3, 5, 10, 30, 100 Myr converted into the ACS photometric 
system using the solar metallicity atmosphere models from the AMES library (solid line), the NextGen models of solar metallicity (dashed 
lines) and the same of LMC metallicity Z = 0.01 (dot-dashed lines). Left panel: Representation of the displacements in the CMD of 
identical points on the evolutionary models caused by changing the three grids of atmosphere models used. 

those that fall blue-ward of the minima between the two 
gaussian fits of the CMD of LH 95, shown in Figj4TZe/£). 
This second selection border is drawn on the CMD of 
the area of LH 95 in Fig. [6j we prolong this border at 
lower and higher luminosities, as shown with a dashed 
line in Figure [6] At the faint end a linear extrapolation 
in the CMD is chosen; for the bright end we arbitrar- 
ily define a cut in order to remove a handful of stars 
at I ~ 22 mag, undoubtedly located on the MS. There 
are 162 stars, which we consider to be field subtraction 
remnants, plotted with gray dots to the left of the MS 
border in Fig. [6j These sources will be excluded from 
our subsequent analysis. 

We stress that, despite the improvement of this field 
subtraction refinement with respect to that of Paper I, 
some uncertainty in the membership assignment on the 
tail of the distributions (close to the cut shown in Figure 
[6]) are still present. This effect, which might affect the 
attribution of the oldest stars, should not be significant. 
The population, within the central region, considered MS 
field is about as numerous as that considered part of the 
young system (1904 and 1942 stars respectively). As a 
consequence, the number of young members erroneously 
removed being too blue should be roughly compensated 
by a similar number of MS stars redder than the minima 
between the two Gaussian (Figure [4]) . This guarantees 
that the number of old PMS members of stars is fairly 
accurate, but could lead to an underestimation of the 
ages of some of those. In any case, since in this work we 
demonstrate a significant age spread for LH95, this fact 
would just reinforce our findings. 
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Fig. 3. — The stellar population of LH 95 with, overlaid, the 
FRANEC models converted into observable colors and magnitudes 
by extrapolating the results of the AMES model atmospheres for 
Z = 0.01. 

the control field than in LH 95, in agreement to what 
we discussed above. This difference demonstrates quan- 
titatively that the contamination of the CMD of LH 95 
by field stars cannot be 100% removed, and this applies 
mostly to the reddest part of the field MS in the area of 
LH 95 due to higher reddening. As a consequence we fur- 
ther remove from the sample of cluster-member sources 



5. AGE DETERMINATION OF PMS STARS 

In this section we apply our new maximum-likelihood 
technique for the derivation of the ages of the PMS stars 
of LH 95 from their photometry. This method allows the 
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Fig. 4. — Left: Hess diagram of the stellar density in the area of LH 95, with no field subtraction. There are two distinct stellar sequences 
apparent, the MS attributed to the LMC field and the native PMS population. The results of the simultaneous fit of two Gaussians along 
diagonal cross-sections are drawn, using filled circles to indicate the centers of the Gaussians for every strip, open circles for the minima 
between them, and white lines following the directions of the strips for the standard deviations of the Gaussians. Right: The same as in 
the left panel for the area of the control field. This CMD demonstrates the complete lack of any PMS population in the local LMC field, 
making the distinction of the true PMS stellar-members of LH 95 quite straightforward. A single Gaussian is used to fit the MS of this 
CMD along identical cross-sections as for LH 95. 
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Fig. 5. — CMD broadening of the MS in color as a function of 
luminosity for the area of Lfl 95 (solid line) and the control field 
(dashed line). The smaller broadening of the MS in the control field 
than in LH 95 demonstrates that field stars could not be completely 
removed within the application of the statistical field-subtraction 
in Paper I. Here we refine our selection of true stellar-members for 
LH 95, based on the results of Gaussian fits to the distribution of 
the MS in the field and the system. 

determination of the age of a PMS star from its position 
in the CMD in a probabilistic manner, accounting for the 
physical and observational effects that affect the colors 
and magnitudes of all observed PMS stars. 

In order to estimate the age of a young stellar system 
when only photometry is available, it is common practice 
to simply overlay a set of isochrones to the CMD and, by 
interpolation, assign an age to all the members. This ap- 
proach, howeve r, can lead to inaccurate results ( Goulier- 
mis et al. 2010). This is because every star is scattered 
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Fig. 6. — Refinement of the removal of the field contamination 
in the CMD of LH 95. The drawn solid line follows the minima of 
the best-fitting Gaussians to the MS and PMS stellar distributions 
before any field subtraction (see Fig.Hje/i), showing that the bluest 
sources of the population clearly belong to the MS. The dashed 
line is our chosen continuation of the previous, at higher and lower 
luminosities. Sources located to the blue of this MS border and 
its extensions (dashed lines) are excluded from our analysis as true 
field-members. 

from its original position in the CMD by a number of fac- 
tors, and specifically by photometric errors, differential 
reddening, unresolved binarity, stellar variability, non- 
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photospheric excess due to accretion, and scattered light 
from circumstellar material. 

These effects, when not accounted for, bias the ob- 
served age distribution, both increasing the observed 
spread and systematically shifting the mean cluster age. 
These factors can be disentangled only on the basis of 
individual stars, which is the case when spectroscopic 
measurements are available for all stars. However, such 
observations on large numbers of stars are not technically 
feasible, and photometric observations are the most ap- 
propriate for collecting data on vast stellar numbers, like 
in the case of LH 95. In this case, a convenient way to 
estimate the overall effect of observational and physical 
biases to the original CMD positions of PMS stars is to 
compare the observed data with synthetic CMDs, which 
are constructed by the application of these biases to one 
or more theoretical PMS isochrones. There are some ex- 
amples of such an approach presented in the literature, 
especially for the investigation of the relation between 
the observed luminosity spread and the true age-spread 
among PMS star-members of young stellar systems (e.g., 
Burningham et al.||20051 1 Jeffries et al.||2007 |Hennekem-| 
per et al.||2008NHillenbrand et al.||2008[). 



TABLE 1 

IMF SLOPES CORRECTED FOR BINARITY, 
AS A FUNCTION OF BINARY FRACTION 



In this context, |JNaylor &; Jettries| ( 2006 1 (hereafter 
NJ06) presented a rigorous approach to derive the age 
of young stellar systems accounting for photometric er- 
rors and binarity. They present a maximum-likelihood 
method for fitting two-dimensional model distributions 
to observed stellar CMDs. In their application, the two- 
dimensional models are probability density distributions 
in the CMD, obtained by applying binarity to a PMS 
isochrone of a given age, considering that the model rep- 
resenting a coeval population in the CMD is not a line - 
like a theoretical isochrone - but a 2D density distribu- 
tion. There are several advantages in the application of 
this approach. First, it enables one to account rigorously 
for the photometric errors in both color and magnitude 
when the fit is performed. Second, binarity is considered 
in the fitting process itself. Finally, this method allows 
the evaluation of its goodness-of-fit. 

In the present study we further develop the method 
of NJ06, by extending its application to a realistic case 
of a PMS population, including in the computation of 
the two-dimensional models not only binarity, but also 
variability, crowding, and differential extinction. 

5.1. From isochrones to density distributions in the 

CMD 

As in NJ06, we wish to construct a stellar density dis- 
tribution in the CMD for every age in the considered set 
of evolutionary models. This distribution represents the 
probability for a star of an assumed age to be found at 
a given point in the CMD. We consider the two families 
of evolutionary models (FRANEC and Siess) that we 
converted int o ma gnitudes in the HST/ACS photomet- 
ric system in § |3.1| For every considered age, we populate 
the appropriate isochrone with 5T0 6 simulated stars. We 
assume the IMF of LH 95 derived in Paper I and correct 
it for binarity. The change in the measured IMF slopes 
caused by binarity was already investigated in Paper I, 
through a Monte-Carlo simulation. Here we adopt these 
results, summarized in Table [T] We found that the IMF 
for intermediate and high-mass stars (M > 1 M Q ) is al- 
most insensitive on binaries, whereas for low-mass stars, 
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Note. — In th e units we use t o ex- 
press IMF slopes, a| Salpeter||l955l l IMF 
has a slope x = 1.35"! 

t Measured IMF from Paper I. 

the slope changes of up to 0.4 units adding an increasing 
fraction of unresolved binaries. Neglecting binaries leads 
to a shallower measured IMF. 

For smaller masses, not having a direct measurement 
of the IMF in LH 95, we assu me a value x = 0.3 compat- 



ible with the Kroupa ( 2001[ ) IMF. This value was also 
assumed in Paper 1 to compute the binarity-corrected 
IMF slopes reported in Table [Tl It is worth noting that 
small changes in the slope of the assumed IMF do not 
influence significantly the age s fo r single stars derived 
using our method (see Section [5~2| ) , but could affect the 
evaluation of the goodness-of-fit (Section [6| we use to 
determine the age spread. 

After simulating a population of stars for every age, 
(for the two families of evolutionary models, and all the 
binary fractions / reported in Table [T] we displace with 
a random process the stars in the CMD according to the 
dictations of several factors. They are discussed below. 

(i) Rotational variability: We simulate statistically the 
variation of magnitudes in the F555W and F814W 
bands due to variability cau sed by stellar r o tation 



We consider the results of Herbst et al. (2002). 
These authors surveyed the PMS population of the 
Orion Nebula Cluster (ONC) in the /-band over 45 
nights, obtaining the light curves for ~ 1500 candi- 
date members. We consider their catalog of peak- 
to-peak variations. A single imaging of a young 
stellar system, like LH 95, statistically measures 
the flux of variable stars in random phases of their 
variation. Considering, at this point, variability 
only due to dark spots on the stellar surface, the 
measured fluxes are equal or fainter than those pre- 
dicted if no variability is considered. Therefore, if 
we assume for simplicity sinusoidal light curves, the 
statistical distribution of /-band variations for a 
star with a measure peak-to-peak A/ max is the dis- 
tribution of A/ max • (sin a; + l)/2 with < x < 2tt. 
We compute this fo r every star of the catalog of 
Herbst et al. ( |2002 ) and combine the derived dis- 
tributions. The result, shown in Fig. [7^, is the 
observed /-band variation distribution that we as- 
sume as valid for a single epoch photometry in a 
young stellar system. The peak at very low A/ is 
mainly due to the large fraction of stars not show- 
ing significant rotational variability. We displace 
the brightness in the F814W filter of each of the 
simulated stars along each isochrone by a value ran- 
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domly drawn from this distribution. 

The amplitude of the V— band variation due to 
dark spots is generally larger than in /, and the 
ratio between the two depends, in general, on the 
temperature of the star-spot, on the temperature of 
the photosphere, and the fill ing factor of the spot . 
According to the results of jHerbst et al.| (fl994), 
the ratio a = AI/AV spans from 0.39 toOTSS. 
The variations in the two bands should be corre- 
lated if the observations were close to co-temporal. 
In our photometry, however, as mentioned in Sec- 
tion [2j V and I observations were taken s eparately, 
one 4 day s after the other. According to Herbst et 
al.| ( |2002[ ), the rotational period distribution of T- 
Taun stars is bimodal, with a peak at P ~ 2 days 
and another at P ~ 8 days. Even considering the 
longest periods, the time delay of our observations 
covers about half period. As a consequence, we 
cannot consider the V— and /—band variations to 
be correlated. Wc therefore simulated the V— band 
variation distribution from the /—band one, mul- 
tiplying each of the simulated / values for a fac- 
tor I /a, where a was randomly chosen within the 
range 0.39-0.88. Then these V variations were ap- 
plied randomly (regardless Al) to all the simulated 
stars. 

(ii) Accretion: On-going accretion from the circumstel- 
lar disk onto the stars also affects the observed 
fluxes. In the optical wavelength range, two com- 
ponents contribute to the accretion excess: an opti- 
cally thin emission generated in the infalling flow, 
and an optically thick emissi on from the heated 
photo sphere below the shock (Calvet & Gullbring 
1998). The optically thin emission is character 



ized by a line-dominated spectrum, similar to that 
of a HII region; in particular the bright Balmer 
emission is commonly us ed to derive mass accre tion 
rates for PMS stars (e.g. JDe Marchi et al.|2010[ ). In 
V— and /— band the optically thick emission domi- 
nates, producing an excess, in spectroscopy reffered 
to as "veiling" , analogous to the presence of a hot 
spot on the stellar photosphere. For a given accret- 
ing star, the spectrum and intensity of the excesses 
may vary with time, due to both variations in the 
mass accretion, and stellar rotation. 

In order to model the overall effect of accretion 
on a large sample of low - mass stars, we utilize the 



results of Da Rio et al. (2010). In this work the 



stellar population of the Orion Nebula Cluster was 
studied from multi-band optical photometry and 
spectroscopy. In particular, for every star, the veil- 
ing excesses were determined from the comparison 
of the observed colors with the intrinsic ones, after 
modeling the color displacements caused by accre- 
tion and extinction. For the ONC, the authors find 
about 60% of the members showing evidences of ac- 
cretion veiling. We consider this stellar sample, and 
restrict to the mass range 0.3M Q < M < 3M Q , rel- 
evant to our study. Then, we displace each of the 
simulated stars for every isochrone of an amount 
AV and AI equal to the veiling of a randoml y 
chosen Orion member from Da Rio et al. (2010). 



An example of the distribution of displacements in 
color and magnitude is shown in Figure [7]3. 

We stress that, statistically, the effect of variabil- 
ity in the mass accretion is al ready included in the 
results of Da Rio et al. (20101; this is because their 
observation produced a "snapshot" of the veiling 
excesses for ~ 1000 members at a given time. It 
could be argued that the members of LH 95, with 
an average age about twice as older as Orion, could 
be affected by less veiling than for the Galactic clus- 
ter, since mass accretion de creases with stellar age 



( |Sicilia-Aguilar et al.||2006"| . With no direct mea- 
surements of accretion in LH 95 we are unable to 
accurately constrain the real distribution of veiling 
excesses for our region. Nevertheless , there are ev- 
idences (e.g., Romaniello et al.]|2004 De Marchi et 
al. 2010) suggesting that accretion lasts longer in 
the Magellanic Clouds, because of the lower metal- 
licity. Because of this we are confident that the 
distribution of V and / excesses we assumed could 
be appropriate for LH 95. Future measurements of 
mass accretion in our region (e.g., from Ha pho- 
tometry) could help to better constrain this param- 
eter. 



(Hi) Binarity: For every isochrone, we associate pairs 
of stars from the simulated sample. This is done 
randomly, but allowing only configurations within 
a mass range, defined by constraining th e mass ra- 

tio to be a flat distribution from 0.2 to 1 ( Woitas et 



al 



2001), and do not allow binary systems to form 
h mass ratio < 0.2. We let the binary fraction be 
a free parameter, and compute different sets of two- 
dimensional models for different choices of its value. 
This is because, as we discuss later in detail, the 
assumed binary fraction influences significantly the 
derived ages, and therefore we present and compare 
findings obtained with different choices of this pa- 
rameter. We assume a unique binary fraction and 
mass ratio distribution for the whole PMS mass- 
range. In principle high-mass stars are known to 
have both higher binary fraction and higher binary 
mass ratio than intermediate and low-ma ss stars 
(e.g., |Lada||2006| |Zinnecker & Yorke||2007|) . How 



Lada 


2006 


this c 


oes i 



our analysis of stellar ages we only consider PMS 
stars and not the 'evolved' high-mass MS stars of 
LH 95. 

(iv) Differential extinction: In Paper I we derived the 
extinction distribution for LH 95, measuring the 
reddening of the upper main-sequence (UMS) from 
the ZAMS. In that work not all the UMS were con- 
sidered to derive the Ay distribution, but some of 
the reddest were removed, since we suspected them 
to be Herbig Ae/Be stars. The average extinction 
is quite low, equal to A555 ~ 0.6 mag, and the its 
distribution spans from A 555 = to A 555 ~ 1. We 
use the relation ^4555/^814 = 1-85, as it was com- 
puted in Paper I for the AC S photometric s ystem , 
using the extinction law of Cardelli et al. (1989). 
In our simulations we apply a random value of ex- 
tinction to each star, drawn from the real measured 
A555 distribution of LH 95. 
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Fig. 7. — Upper panel: Distribution of the expected mgi4— band 
variation due to dark spots and stellar rotation in single e poch pho- 
tometry, o btained from the peak-to-peak distribution of |Herbst et| 
|al.| ( |2002} . Lower panel: Displacements in color and magnitude 
caused by accretion veiling, from the findings of |Da Rio et al.| 
( [2010] ). 

(v) Source confusion: LH 95 is a loose association, not 
centrally concentrated, which comprises at least 3 
subclusters of somewhat higher stellar densities. 
As a consequence, source confusion due to crowing 
does not play a dominant role. Nevertheless, we in- 
clude in our simulations a 5% source confusion, by 
adding the flux of this fraction to sources regardless 
of their assumed binarity, i.e., by allowing for ran- 
dom confusion also between binaries and between 
binaries and isolated stars. We stress that, since 
this additional pairing is done randomly because 
of the IMF the modeled source confusion is also 
magnitude dependant, favoring the disappearance 
of faint systems. 

For the computation of the 2D stellar densities we do 
not include photometric errors as additional sources of 
spread. While the aforementioned biases, which we in- 
clude in the simulations, can be considered as part of 
the model, the photometric errors are part of the data, 
which the model is fitted to, and therefore photometric 
u ncert ainty will be consider later in the fitting process 
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Fshould be finally noted that the constructed 2D dis- 
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Fig. 8. — Example of a 2 Myr isochrone from the Siess mod- 
els, which is converted into a 2D density CMD distribution (blue 
contour) , after the application of binarity (for an assumed binary 
fraction / = 0.5), accretion , var iability, differential extinction and 
crowding, as described in § |5.1| The original isochrone is overlaid 
with the red line. The dots are the stars of our photometry with 
error bars indicating the lcr photometric errors. An average red- 
dening of E(V — I) = 0.275 mag, typical for LH 95 (Paper I), is 
applied. 

tributions are not re- normalized to, e.g., a total inte- 
gral of 1, because isochrones of different ages have dif- 
ferent maximum mass. Instead, we use the same sam- 
ple of 5 million simulated stars with masses drawn from 
the assumed IMF and we consider only the fraction in- 
cluded within the limits of each isochrone. In this man- 
ner, although the total number of synthetic stars changes 
slightly among different 2D distributions, the number of 
stars per mass interval is constant. In other words, since 
our method to assign stellar ages and determine the clus- 
ter age is based on relative comparison of probabilities, 
an overall normalization of all the models is completely 
unnecessary. 

In Fig. [8] we show an example of the result of the appli- 
cation of our simulations to a 2 Myr isochrone from the 
Siess family of models. In this figure it is evident that the 
CMD-spread of a 'perfect' coeval PMS stellar population, 
produced by the biasses discussed above, can generate an 
observed dislocation of the stars of more than I mag in 
luminosity and up to I mag in color. The bifurcation 
seen at the faint end of the simulated distribution is an 
effect of binarity. This is due to the fact that the evolu- 
tionary models do not cover masses lower than 0.1 M©. 
As a consequence our random selection of binary pairs is 
unable to simulate binaries with the second component 
below this limit. This effect does not alter the results of 
our analysis, because, as it is seen in Fig.JHJ our photom- 
etry does not reach these masses. 
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function of the age t as: 



.5 7.0 
og (age) 
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Fig. 9. — Example of likelihood functions, for 3 stars from 

our photometric sample. The maxima of the functions indicate the 
assigned best-fit ages for the sources, equal to 1.3, 5.5 and 19 Myr 
respectively in this example. 



5.2. Derivation of stellar ages. The single star case. 

In this section we use statistics to address the issue 
of deriving the most likely age of a single PMS star, for 
which the only information is its measured brightness 
and the corresponding photometric errors. Specifically, 
we derive the age of each star in our catalog using a 
maximum-likelihood method by comparing the observed 
position of each star in the CMD with a complete set 
of simulated 2D stellar distributions, constructed from 
theoretical isochrones with the method described above. 
Each of these simulated distributions provides an esti- 
mate of the probability of finding a star in a specific 
CM D lo cation, providing that we know its age. In Sec- 
tion |5.4| we will then generalize this approach to derive 
the best age of the entire cluster, considering all the mea- 
sured stars simultaneously. 

In practice, if the photometric errors were zero, the 
most likely age would be that of the modeled distribution 
that gives the highest probability in the CMD position of 
the star. In reality, since there are errors in the measured 
?7i555 and Ti8i4 magnitudes of the star, the observed po- 
sition of the star in the CMD is also probabilistic, defined 
by the gaussian distributions of the photometric errors. 
Therefore, if we consider a star observed in the position 
(Vi,Ii) on the CMD and a model stellar distribution of 
a given age t, for every point (Vj,Ik) of the CMD-space 
there are two probabilities that should be specified: a) 
the probability, defined by the 2D model, that a star of 
age t is located in (Vj,Jk), and b) the probability, de- 
fined by the Gaussian distributions of the photometric 
errors, that the star we observe in (Vj, Ij) is instead lo- 
cated in (Vj, Ik) but displaced by the photometric errors. 
The product of these two probabilities, integrated over 
the entire CMD-space, provides the likelihood that the 
considered star has an age equal to t. This likelihood, 
normalized on all the possible ages, represents the exact 
probability that the star has an age between t and t + dt. 

Since the photometric errors in the color-magnitude 
space are correlated, we apply this method in the 
magnitude-magnitude plane. Mathematically, we define 
the likelihood Li(t) for the i—th star in our sample as 



Li(t) 



p(V,I,t) Ui(V-Vi,I-Ii)dVdI (1) 



v,i 



where p{V, J, i) is the surface density computed from the 
isochrone of age t, and 



Ui{V-V h I-h) = 



1 



27TCT Vi °~ Ii 



■[(- 



(2) 



is the 2D Gaussian representing the photometric error of 
the i—th star. 

In Fig. [9] we show an example of the likelihood func- 
tion Li(t) computed numerically for three stars in the 
LH 95 catalog, assuming 2D models obtained from the 
FRANEC grid and assuming a binary fraction / = 0.5. 
The age assignment for every star is performed by max- 
imizing Li{f). For example, the three stars whose Lj(t) 
functions are shown in Figure |9j will have most probable 
ages of 1.3 Myr, 5.5 Myr and 19 Myr. 

In order to estimate the accuracy in the age determina- 
tions, we isolate for every likelihood function the range, 
in log (i) that includes 68% of the area under the curve, in 
analogy with the ±lcr interval for a normal distribution. 
Since for some stars the shape of Li (t) is highly skewed, 
we choose to define the limits of the interval to delimit 
identical areas of 34% of the total to the left and to the 
right from the peak of the distribution. We then call e 
the width of the range, in units of log (t) . Indicatively, if 
a Li (t) we re normal, o~(t) would be equal to e/2 

Fig. Mleft) 
in the (JMD, 
FRANEC and 

Fig. |10[ riff/7i) we show the corresponding variation of 
e along the CMD. It is evident that, whereas the age 
assignment turns out to be relatively precise in the 
PMS regime, for higher masses, close to the PMS-MS 
transition (msu < 21) the uncertainty increases up to 
e ~ 1.5 dcx. This is due to the fact that PMS isochrones 
of different ages are closer to each other in this range, and 
tend to overlap in the MS regime. The higher values of e 
on the blue side of the sequence are also due to the fact 
that theoretical isochrones of ages > 10 are also close to 
each other, since PMS evolution is significantly slower for 
these ages than for younger stages of PMS contraction. 

We compare the distributions of the best-fit ages for 
our photometric sample changing the binary fraction and 
the original family of evolutionary models. This is shown 
in Fig. [TT] where we also indicate the average age for each 
distribution. It this figure is evident that increasing the 
binary fraction leads to the pred i ction of older ages. This 



shows the derived best-fit stellar ages 
assuming the evolutionary models of 
a binary fraction / = 0.5. In 



confirms the findings of Naylor 



(2009). Indeed, a higher 



fraction of binary systems produces a simulated 2D dis- 
tribution shifted towards higher luminosities, hence older 
measured ages. Passing from / = to / = 0.8, we pre- 
dict a difference in age of up to 0.2 dex (which is 60% 
older ages.). 

In the same figure we also overlay the distribution of 
ages computed only by interpolation between theoretical 
isochrones, without the addition of any observational and 
physical bias except for the average extinction of LH 95. 
In this case we find that the average age is higher than 
that derived from our statistical method assuming / = 0. 
This is due to the effect of dark spots and accretion. The 
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Fig. 10. — Left: Color-magnitude diagram with the PMS stars color-coded according to their individual ages derived with our maximum- 
likelihood method. The model distributions were derived with the use of theoretical isochrones from the FRANEC family of mo dels. The 
assumed binary fraction is / = 0.5. Right: The same CMD with the stars color-coded according to their e parameter (see § |5.2| 

first biases the predicted luminosities towards fainter val- 
ues; the second towards brighter luminosities and bluer 
colors. In particular, the color effect caused by veiling 
affects more the position in the CMD than the increase 
in the fluxes, requiring an older isochrone to fit the ob- 
served sequence. 

In conclusion, we find that, in the derivation of stellar 
ages for a young PMS cluster, neglecting binaries leads 
to an underestimation of the ages, neglecting variability 
and accretion to an overestimation of those. The two 
effect counterbalance, in our case, for a binary fraction 
/ = 0.5 and assuming star-spot and accretion values from 
the 2 Myr Orion Nebula Cluster. 
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5.3. Spatial variability of stellar ages 

As we mentioned in Section [2j LH 95 is not a centrally 
concentrated cluster. It presents 3 main concentrations 
of PMS stars, named in our previous work "subcluster 
A, B, C", located respectively at the eastern, western, 
and southern edges of LH95. They include about half 
of the stellar population of the region, while the rest of 
the population is distributed over the area between the 
three. Since the projected distance between the stellar 
clumps is quite large (15 - 20 pc), it could be plausible 
that star formation did not occur simultaneously all over 
the entire region, but sequentially. A Galactic ex ample of 
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such scenario can be found in the Orion complex ( Briceno 
2008|), where different populations, from 2 to 10 Myr old, 



are distributed over a similar range of distances. 

We investigate if the stellar sub-groups in LH95 present 
different average ages, in order to rule out the possibility 
that the meas ured broad age distribution for the entire 
region (Figure 11) is simply due to a superposition of 
populations witn different ages, instead of the effect of 
observational biases and, maybe, a real age spread. To 
this purpose we isolate the members in the three sub- 
clusters, using the same selections as in Paper I. Then 
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Fig. 11. — Distribution of maximum-likelihood stellar ages for 
the LH 95 members, computed using 2D CMD distribution mod- 
els, based on FRANEC (upper panel) and Siess (lower panel) evo- 
lutionary models. We assume different binary fraction values / 
from 0.0 to 0.8. For each family of isochrones, we also show the 
distribution of ages computed only by interpolation of isochrones 
in the CMD, without accounting for physical and observational bi- 
ases (black lines) . The vertical lines (colored appropriately) in 
each panel indicate the mean value of the age for each plotted 
distribution. 
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Fig. 12. — Comparison between the age distribution of the entire 
region and the three denser sub-clusters, for FRANEC evolution- 
ary models and / = 0.5. 

we consider the distribution of ages limited to each of 
the three subregions, in comparison wit h th at of the en- 
tire LH 95. This is illustrated in Figure [12] It is evident 
that, in each region, the distributions are very similar, 
spanning from less than 1 Myr to 25 Myr, with a peak 
at ^5 Myr. The only evident anomaly is an overabun- 
dance of young members for sub-cluster B. However, this 
group of stars encloses the densest concentration of PMS 
stars in the whole region, with more than 150 star within 
~ 2 pc. Because of this, stellar confusion due to crowding 
should be particularly dominant in this area, producing 
a number of overluminous stars due to unresolved su- 
perposition, and, thus, a fraction of apparently younger 
members. This is, however, not a problem for the analy- 
sis of the clu ster age (Section |5.4[ ). In our CMD modeling 
(Section [OJ we have purposefully included the effect of 
crowding affecting 5% of members. The overabundance 
of young stars in subcluster B is caused by about 40 stars, 
which are ~ 2% of the total stellar sample. 

Considering this, we exclude the presence of any spatial 
variations of the ages in LH 95. Despite the irregular 
density of stars and the relatively large size of the region, 
LH95 must have formed as a whole, in a unique episode 
of star formation. 

5.4. Age of the entire cluster 

In section [572] we derived the most likely age for every 
star in the sample. Here we derive the best-fit age for 
the entire stellar system. In a similar way as in NJ06, we 
compute this age by maximizing the "global" likelihood 
function L(t) — J\iLi(t). To compute L(t) we restrict 
our treatment to the stellar sample with 22 < TO555 < 27. 
In particular, we want to eliminate any biases introduced 
by detection incompleteness, which would lead to the de- 
tection at the faint-end of more young bright members 



than old faint ones. The detection completeness reaches 
~ 90% at rri555 = 27 (Paper I), and we consider this 
threshold as a rea sonable li mit fo r our analysis here. We 
stress that, as in Nay lor ([2009), we could include the 
detection completeness tunction in the 2D models de 



rived in Section [5. 1| and then recover the implied cluster 
properties using our maximum likelihood procedure con- 
sidering all the stars. In this way, we could test for bias 
also as a function of our chosen completeness limit. Nev- 
ertheless, considering the large stellar sample we have at 
disposal, we prefer to limit out analysis on the luminosity 
range where we that detection incompleteness does not 
affect our photometry. In this way we remove any addi- 
tional uncertainty in the cluster age caused by a possible 
inaccuracy of the measured completeness function. 

In addition, we do not include in the analysis the bright 
members of both the MS and the PMS-MS transition, 
for two reasons: 1) as shown in Fig. [To| right), the ac- 
curacy of the age estimation for these sources is poor, 
and 2) these sources have very small photometric errors, 
which means that the 2D Gaussian distribution repre- 
senting the errors (Eq. [2l goes to zero very fast in the 
neighborhood of the specific CMD loci, and so does Li(t). 
Such very small values of Li (t) would dominate the prod- 
uct LJ Li (t) , biasing the resulting global likelihood func- 
tion. From a mathematical point of view this effect is 
legitimate, because it only means that during the fit the 
data points with the smallest uncertainties dominate the 
statistics. However, in practice if the observational er- 
rors are much smaller than the intrinsic uncertainty of 
the models (in our case, the uncertainty in the theoret- 
ical evolutionary model computation and the modeling 
of the scattering) a strict application of the former can 
lead to unrealistic results. 



In Fig. 13 wc show the derived L(t) functions, in log- 
arithmic units, for both FRANEC and Siess isochrones 
and different choices of binary fraction /. The plateau 
at old ages, in comparison to the cut-off at young ages, 
is due to the fact that old isochrones tend to be closer 
to each other and closer to the sequence of PMS stars 
observed in the CMD than the young isochrones. As a 
consequence, the 'distance' from a given star, in units of 
er of the photometric errors, is on average lower for the 
older isochrones than for the younger, resulting to higher 
probabilities for the former. 

Results of our analysis are given in Tabled where wo 
show the ages at the peak of each L(t), as wellas the asso- 
ciated values of likelihood for each choice of parameters. 
Our analysis also shows that the best-fit age is very sen- 
sitive to the assumed bi nary fraction, consistently with 
what we found in Section 5.2 analyzing the distribution of 



ages for single stars, as it increases of up to 60% with in- 
creasing binarity (depending on the considered isochrone 
models). The peak values of logL(i) given in Tabic [2] do 
not have an absolute meaning, given the arbitrariness of 
the normalization of the 2D density models p(V,I,t). 

In Fig. [14] we show again the observed CMD, with 
the best-fit 2D density models from FRANEC (left) and 
Siess (right) grids of isochrones superimposed for / = 0.5, 
for the luminosity range we considered in our analysis. 
The two models are very similar, with FRANEC pro- 
ducing a slightly over-density at the reddest end of the 
considered magnitude range. 

From Figure[14j it is evident that an additional spread, 
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Fig. 13. — Global likelihood functions L(t), derived from the 
simultaneous fit of the the 2D stellar distribution models to all 
PMS stars in our photometric sample with 22 < 771555 < 27 based 
on FRANEC (upper panel) and Siess (lower panel) evolutionary 
models. We assume different binary fraction values / from 0.0 
to 0.8. The vertical lines (colored appropriately) in each panel 
indicate the age of the maxima of each curve. 



TABLE 2 

Most-likely ages for the entire cluster. 





FRANEC Models 


Siess Models 


/ 


age (Myr) 


log L max 


age (Myr) 


log L max 


0.0 


2.8 


2615.26 


3.0 


2492.16 


0.2 


3.0 


2578.42 


3.3 


2456.14 


0.4 


3.6 


2548.63 


3.6 


2431.42 


0.5 


.3.9 


2533.37 


3.8 


2415.62 


0.6 


i.l 


2517.49 


4.0 


2398.78 


0.8 


4.4 


2492.71 


4.3 


2377.05 



besides that modeled in the 2D density distribution, is 
evident in our data. In particular, there are observed 
members at the sides of the sequence which may not be 
well modeled by the simulated distributions for a sin- 
gle cluster age. This may be indicating of an additional 
age spread in LH95. In the next sections we study this 
possibility. 

6. AGE-SPREAD OF PMS STARS IN LH 95 

6.1. Verification of an age-spread 

In the previous section we have determined the aver- 
age age of LH 95 isolating the 2D model that provides 
the best fit to the observed sequence, for two families of 
theoretical evolutionary models. In this section we wish 
to investigate whether the spread of the PMS stars in 
the observed CMD is wider or not than that of the sim- 
ulated 2D density distributions as they are produced by 



acc ount ing for the sources of CMD-broadening discussed 
in § 5.1 after we additionally apply the photometric un- 
certainties. If indeed the observed CMD broadening of 
PMS stars is wider than the synthetic derived from a 
single isochrone, under the hypothesis that we do not 
have underestimated the sources of apparent luminosity 
spread in the CMD, this would suggests the presence of 
a real age-spread in LH 95. The method applied in this 
section, described below, also allows us to qualify the 
goodness of our fit. 

This method can be compared to a standard y 2 min- 
imization, in which the best fit model is obtained by 
minimizing the \ 2 , and this value is then compared to 
the expected distribution of \ 2 given the degrees of free- 
dom. A measured \ 2 which significantly exceeds the ex- 
pected distribution indicates a poor fit, while a value too 
low suggests that the fit is suspiciously good, which is 
the case, e.g., of overestimated errors. In our statistical 
framework, the procedure is qualitatively similar, with 
two main differences: i) since the measured log L(t) is 
maximized, a poor fit is found when the measured value 
is lower than expected; ii) the expected distribution of 
logi(tbcst-fit) is not analytical but can be derived nu- 
merically. In general, this distribution varies for different 
assumed 2D models, and depends on the number of data 
points. 

We compute the expected distribution of log L as fol- 
lows. We consider the two best-fit 2D models shown in 
Fig. |14[ limiting ourselves, for the moment, to the case 
of an assumed binary fraction / = 0.5. We limit the 
model to the same range used to derive the cluster age 
(22 < m555 < 27), and simulate a population of stars 
randomly drawn from the distribution of the model itself, 
as numerous as the LH 95 members in the same luminos- 
ity range. We consider each of the simulated stars, and 
displace its position in the CMD by applying the average 
photometric error measured in its immediate CMD area. 
In particular, for each simulated source we consider the 
60 observed stars closer to it in the magnitude-magnitude 
plane. Than we randomly assign the photometric error 
of one of these 60, with a weight proportional to the in- 
verse of the distance in magnitudes. We then apply the 
same method described in § |5.4| to the simulated popula- 
tion, deriving log L. We iterate this procedure 100 times, 
each time randomly simulating a different population of 
test-stars from the same best-fit 2D model. In this way 
we obtain 100 values of log L describing the expected dis- 
tribution of this variable for the considered model, which 
is compared to the measured value computed for the ob- 
served sequence of PMS stars. 

For the two models shown in Fig. [Mj i.e., the 
FRANEC 3.9 Myr and the Siess 3.8 Myr models for 
/ = 0.5, we obtain predicted logL values of 2739 ± 20 
and 2473 ± 19 respectively (±lcr interval). In compari- 
son, the measured values of 2533 and 2416 (see Table [2| 
are significantly lower, at ~ — IOct and — 3a from the ex- 
pected values. This implies that the observed sequence 
cannot originate from the bcst-Rt models by additionally 
applying photometric errors, but that the stars tend to 
be statistically located in regions of the CMD where the 
2D model has lower density p (see Eq. [T]) , and hence a 
lower average log Li for the stars. For guidance, the low 
density regions of the 2D models in Fig. [Ml are those 
indicated with lighter blue colors, i.e., at lirgh luminos- 
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Fig. 14. — The observed CMD with the 2D density models that provide the best fit to the observed sequence of PMS stars, assuming 
binarity of / = 0.5 for the FRANEC (left panel, age 3.9 Myr) and the Siess (right panel, age 3.8 Myr) family of evolutionary models. 
The dashed lines denote the subsample of the photometry used in this analysis. For clarity, only the portion of the 2D models within this 
luminosity range is shown. 

ity, where the IMF predicts fewer members, and at both 
sides along the 2D distribution in the CMD, at large dis- 
tance from the peaks of the 2D sequence itself. Since 
the IMF we use to compute the 2D models has been 
derived from the actual observed data, we exclude the 
possibility that the low measured value of logL is due 
to an overabundance of intermediate mass stars in com- 
parison to low-mass stars. As a consequence, our results 
suggest that the observed sequence of PMS stars in the 
CMD is broader than that predicted by the 2D model for 
a single age. There are two possible interpretations for 
this: a) LH 95 is formed by a coeval population, and the 
extra spread in the observed CMD is due to additional 
sources of apparent broadening in the CMD besides those 
we have included in our models; b) LH 95 hosts a real 
age-spread among its PMS stars. From hereafter we will 
assume the second possibility, and, under the validity of 
our modeling, we determine the age spread of the popu- 
lation. 
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6.2. Evaluation of the age-spread 

In order to quantify the age-spread we include in our 
construction of the 2D model distributions, apart from 
the sources of CMD-broadening of PMS stars previously 
considered (§ |5.1[ ), also a distribution in ages. Specifi- 
cally, we consider the best fit age previously measured 
for the cluster (Table [2 ) along with a Gaussian distri- 
bution of ages peaked at the best fit age with standard 



deviations varying from <7 agc = O.lMyr to <r a 



3 Myr. 



We construct, thus, broader 2D model distributions by 
essentially ad ding together different single- age models, 
computed in § |5.1[ each weighted according to the Gaus- 



FlG. 15. — Comparison between the measured value of log- 
likelihood and the predicted one, in units of standard deviation 
of the latter, as a function of the age spread. 

sian distribution. For each value of o~ age we compute the 
log L from our data on the newly constructed 2D model, 
as well as the expected distribution of log L, from sim- 
ulated stars drawn from the model, as described above. 
We stress that the expected distribution of log L is not 
unique, but varies for each assumed value of <7 age , being 
derived from the 2D model computed for the exact value 

Of (Tago 

For every value of <7 agc we compare the "observed" 
logi (i.e. derived from the data) with the predicted 
one. As mentioned earlier, if the first is significantly 
smaller than the second, it indicates a poor fit; vicev- 
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Fig. 16. — The observed CMD with the 2D density models that provide the best fit to the observed sequence of PMS stars, including 
an intrinsic age-spread, for an assumed binarity of / = 0.5. Left: The model obtained from FRANEC isochrones, with an average age of 
3.8 Myr and a Gaussian age-spread of <r agc = 0.18 Myr. Right: The model obtained from Sicss isochrones, with an average age of 3.9 Myr 



and <r a 



0.12 Myr. The dashed lines denote again the subsample of the photometry used in the analysis. 



ersa, a significantly larger value implies a suspiciously 
good fit, analogous to a \ 2 "C 1. The standard deviation 
of the predicted log L distribution allows us to quantify 
these differences in a meaningful way. This is simplified 
by the fact the predicted logL is normally distributed, 
due to the central limit theorem, since log L is the sum 
over a large sample of stars of the individual single-star 
likelihoods log Lj . 

In Figure [15] we present the difference between the 
observed and predicted log L, normalized on the width 
of the predicted distribution, as a function of the as- 
sumed age spread er agc , for the two families of evolution- 
ary models and assuming / = 0.5. The values of <7 age 
for which the difference in the y-axis is zero provide the 
best-fit value of the age spread, the and the projection 
on the x-axis of the ±ler (dotted lines in the Figure), 
the associated uncertainty. We find that, for FRANEC 
models the measured age spread is cr agc = 1.8 Myr 
(FWHM~ 4.2 Myr), while for Siess models this decreases 
to cr age = 1.2 Myr (FWHM~ 2.8 Myr). In both cases the 
uncertainty in cr ago is about 0.2 Myr. 

The observed CMD, with the best-fit 2D density mod- 
els from FRANEC (left) and Siess (right) grids of 
isoc hron es superimposed is shown again for / = 0.5 in 
Fig. 
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but with the best-fit 2D models constructed 
by imduding the measured age-spread of 0.5 dex. In 
this figure is evident how the model density probabilities 
now follow quite nicely the observed broadening of the 
PMS population of LH 95, within the magnitude range 
22 < ?7i555 < 27, which is considered for the evaluation 
of the age and age-spread of the system. The difference 
of these plots with those of Figure 14 where the best-fit 



2D models were constructed assuming coeval PMS pop- 
ulations, is significant. By adding an age spread, the 
model distributions in the CMD now reproduce well the 
observed sequence. 

Our method predicts a 50% higher age spread assuming 
FRANEC models than assuming Siess models, and the 
difference is larger than the uncertainties of each of two 
values. This is a consequence of the particular shape of 
the FRANEC young isochrones for the lowest masses, 
where the modeling of the deuterium burning leads to a 
brighter predicted luminosity. This can be seen already 
in the original isochrones in the H-R diagram (Figure [T]), 
and is re spo nsible for the over-density seen in both Figure 
14 i and \TE^ l at V — I ~ 3 and / ~ 24. In this region of 
trie CMDthe models predict more stars than actually 
observed, regardless the assumed age spread. Since the 
observed population is, therefore, statistically populates 
more the positions of the CMD where the density of the 
2D models is lower, the measured likelihood turns out to 
be smaller than the predicted one even for the "right" 
age spread. This leads to an overestimation of the value 

of CT ago . 

We remark that this is a general behavior our the 
method we have proposed: our technique may overes- 
timate the overall age-spread. This occurs in the case 
where the isochrones do not follow the observed sequence, 
but have a systematically different shape with respect to 
it. In fact, the condition we imposed to constrain the 
measured tr age is based on the requirement that our data 
and a simulated population of test stars obtained from 
the best 2D density model (including an age spread) pro- 
duce similar values of "global" likelihood logL. If the 2D 
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TABLE 3 

Age and age spread for different binary 
fractions 







FRANEC Models 


Siess Models 




/ 


age 


CT agc 


age 


°"agc 






(Myr) 


(Myr) 


(Myr) 


(Myr) 


/ 


= 0.4 


3.6 


1.4 ± 0.2 


3.6 


1.0 ±0.2 


/ 


= 0.5 


3.8 


1.8 ±0.2 


3.9 


1.2 ±0.2 


/ 


= 0.8 


4.1 


2.3 ±0.25 


4.0 


1.4 ± 0.2 



TABLE 4 

Variation of ages and age-spreads with luminosity for 
assumed binarity / = 0.5 





FRANEC Models 


Siess Models 


magnitude 


age 


a age-spread 


age 


a age-spread 


range 


(Myr) 


(Myr) 


(Myr) 


(Myr) 


22 < m 5 55 < 23 


0.9 


2.2 ±0.7 


1.0 


2.8 ±0.7 


23 < m 555 < 24 


1.4 


1.7 ±0.5 


1.8 


2.4 ±0.9 


24 < m 55 5 < 25 


2.6 


0.6 ±0.6 


2.7 


0.9 ±0.8 


25 < m 55 5 < 26 


3.6 


0.8 ±0.4 


3.6 


0.6 ±0.4 


26 < m 55 5 < 27 


1.6 


1.0 ±0.2 


1.7 


0.0 ±0.8 



models have a shape which strongly differs from the ob- 
served sequence, the model would never reproduce the 
data, but there will be always an assumed value of age 
spread which leads to similar values of likelihood as for 
the observed sample. In this case, the derived age spread 
might be higher than the real one. As we discuss in the 
next section, we can quantify any systematic differences 
in the shape of the isochrones with respect to the data, 
and estimate the correct age spread, by analyze the age 
properties of the population is smaller luminosity ranges. 

We investigated how the assumed binary frequency / 
affects the derived age spread. We have assumed / = 
0.4 and / = 0.6 and reapplied the method described 
above. We stress that, although we previously described 
our results for a larger range of assumed binary ratios, in 
the mass range relevant for our study, the value o f / is 
not expected to differ much from 50% (Lada 2006). The 
results are presented in Table |3] We find that a higher 
fraction of binaries leads to a larger derived age spread, 
and, on the contrary, a lower / to a shorter <7 agc . Both 
the results are a consequence of the dependence of the 
average cluster age as a function of /: given that the 
distance between consecutive isochrones decreases with 
age. A given spread in the CMD results in a larger age 
spread if the central isochrone is older, and viceversa. 

6.3. Luminosity- dependence of the age and age-spread 

We want to verify the presence of systematic changes 
of both age and age-spread with luminosity. This effect, 
which can be caused by uncertainties in the theoretical 
evolutionary models, as well as in our modeling of the 
broadening effects in the CMD, can lead to an overes- 
timation of the measured age spread. To demonstrate 
this dependence we divide the photometric sample in the 
considered magnitude range of 22 < m 555 < 27 into five 
sub-ranges 1 mag wide each, and we repeat the method 
described in § |5.4| to derive the best-fit age within each 
subsample. We then apply the method of § [6] to derive 
the corresponding age-spreads, simulating a family of 2D 



models with different assumed <7 ag c around the average 
age measured for each magnitude interval. 

The results are given in Table |1J again for both the 
FRANEC and Siess families of evolutionary models for 
an assumed binary fraction / = 0.5. We find a clear 
variation of both average age and measured age spread 
with luminosity. In particular, for the bright members, 
close to the PMS-MS transition, which correspond to 
masses in the range 1 < M/M Q < 2, the derived ages are 
systematically younger than those for fainter luminosities 
(and lower masses). The differences in age are of the 
order of 4 Myr. Also, the measured age-spread tends 
to be significantly larger than the average one for bright 
stars, and smaller for the lowest masses. In particular, 
for 22 < 771555 < 23 the age-spread derived with our 
method is larger than a = 2.5 Myr, corresponding to a 
FWHM exceeding 5 Myr. Considering the very young 
ages measured in this luminosity interval, and the fact 
that, for our modeling, the minimum age included in our 
simulation is 0.5 Myr, this finding implies that the best- 
fit age distribution for the luminous part of the PMS 
sequence is highly skewed, with a peak at ~ 1 Myr and 
a broad tail towards older ages. On the other hand, for 
the faintest end of our stellar sample, our method finds 
that a more modest age spread is needed to produce a 
simulated distribution of stars as broad as the observed 
sequence. 

The differences in the derived age spreads as a function 
of luminosities may be due to an incorrect modeling of 
the 2D isochrones we have used. In particular, some of 
the broadening effects we included in the model might 
be mass-dependant. For example, variability and accre- 
tion might affect more the brightest members, and less 
the smallest masses, whereas we assumed identical dis- 
tributions of these effects along the entire mass range. 
Considering our results, and averaging the age spreads 
at different luminosities, the average ag e sp read for the 
entire stellar sample, derived in Section [6~2| should be a 
fairly well approximation of the real spread. 

If we trust our models, the age-luminosity dependance 
we measure implies that stars in the low-mass regime 
have been formed some Myrs before the intermediate- 
mass members. However, this is not what is observed in 
young Galactic 1-2 Myr-old associatio ns, where no evi - 
dence of highly biased IMF are found ( Hart mannj 2 3 ) . 
An explanation for this age difference among different 
stellar masses is the assumed binarity. An overestimated 
binary fraction can lead to overestimation of the age for 
the low-mass stars, because the large numbers of unre- 
solved low-mass binaries produce systematically wider 
modeled distributions of low-mass stars in the synthetic 
CMDs. We have investigated this effect by computing 
the average age in the five selected magnitude intervals 
(as those in Table Ej), but with 2D modeled distributions 
constructed assuming no binaries (/ = 0). While in this 
case the mass-age dependance reduces, it does not van- 
ish, with a reduced variation in the derived ages of the 
order of 1 Myr. Naturally, this is an extreme case since 
such a very low binary fraction is not realistic, but unfor- 
tunately the true binary fraction in the considered mass 
ranges cannot be directly resolved at the distance of the 
LMCB 

4 At the distance of the LMC, 1 ACS/WFC pixel corresponds 
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We have investigated more in detail if the mass-age 
correlation we have found can be a biased result to an 
inaccuracy in the modeling of the 2D isochrones. In 
particular, since the effect caused by accretion veiling 
is the most uncertain among those we have considered 
(mostly because of the lack of mass accretion estimates 
for LH 95), we have tried to redo the entire analysis pre- 
sented in this paper excluding accretion from our mod- 
eling. In case of no accretion, we find that the age of 
the cluster is older (of the order of 6 Myr, depending on 
the other parameters), the age spread is larger (several 
Myrs), and still, there is consistent mass-age correlation, 
with the very-low mass stars systematically older than 
the intermediate-mass stars. 

Taking into account the above, and if we assume that 
indeed such a dependence of age with mass should not 
be observed, as in Galactic young clusters, we can at- 
tribute the observed dependence of age with mass to an 
inaccuracy of the PMS models we have utilized. Actu- 
ally, the inconsistency between the ages predicted by the 
FRANEC and Siess isochrones suggests that the PMS 
stellar modeling is still far from completely reliable. For 
a metallicity close to solar, typical of the Galaxy, several 
studies have been presented to constrain the ev olution- 
ary models based on empirical isochrones (e.g., Mayne 
et al. 2007). Certainly new data that sample different 



( |Palla fc Stahler|2000||Palla et al.|20"05{|Tan et al.|2006| ; 

the loose, not centrally concentrated geometry of Lrl 95, 



young associations and clusters in the metal-poor envi- 
ronment of the MCs will be exceptionally useful to test 
the agreement of data with the PMS models also for lower 
metallicities. 

7. DISCUSSION 

We stress that our results on the age determination of 
LH 95 and the evidence for an age-spread in the system 
are valid under the assumption of the correctness of the 
so urce s of spread in the CMD considered in our models 
(§5.1). It should be noted that having included, apart 



from photometric uncertainties, several important bias- 
ing factors, namely differential reddening, stellar vari- 
ability, accretion binarity and crowding, our modeling 
technique is more accurate than similar previous studies, 
which completely neglect these biases. Using our mod- 
els, we cannot reproduce the observed sequence assuming 
a coeval population, but an additional spread in age is 
measured. Nonetheless, as we have shown, our results 
are very sensitive on the assumptions used to model ob- 
served populations in the CMD. In particular, changes in 
the assumed binary fraction alter significantly both the 
average cluster age and the derived age spread. Also, as 
discussed previously, the overall effect of accretion veiling 
influences both this quantities. For our modeling, we as- 
sumed accretion and variability measurements from the 
Orion Nebula Cluster. If on-going accretion in LH 95 is 
weaker than we assumed, for instance because of an ear- 
lier disruption of the circumstellar disks caused by the 
flux of the early-type st ars present in the r egion (3 of 
which are O-type stars, ( Da Rio et al.|[2010c )) or simply 
because of the older age, the measured age would turn 
out be older and the age spread larger that the results 
we have presented. 

Slow star formation scenarios have been proposed 
based on both observational and theoretical evidence 

to ~ 2400 AU. 



and its relatively low stellar density and mass (see Pa- 
per I) wou ld be compatible with a slow conversion of gas 
into stars (Ke nnicutt||1998[ ). Unfortunately, the absence 
of dynamical information on our region, as well as a pre- 
cise estimate of the gas mass do not allow us to infer 
much about the past evolution of LH 95. 

On the other hand, additional unconsidered effects that 
can displace PMS stars in the CMD may still be present 
and affect the derived age-spread; we briefly discuss them 
here. The presence of circumstellar material can cause 
star light to be scattered in the line of sight direction 



(e.g. Kraus fc Hillenbrand 2009 Guarcello et al. 2010) 



This effect may either increase the observed optical flux 
with respect to that of the central star, or reduce it, if 
the central star is partially obscured by the circumstel- 
lar material. This scenario, however, requires both the 
presence of a significant circumstellar environment and a 
particular geometrical orientation, so, if present, would 
affect only a frac tion of the stars in a cluster. In any 



case, as shown in Guarcello et al. (12010b, stars observed 
in scattered light appear significantly bluer than the PMS 
sequence. If such sources are present in our photometric 
catalog, they would affect only the old tail of our age 
distribution, or even be excluded as candidate MS field 
stars. Thus, we are confident that neglecting this effect 
should not affect significantly our findings. Another un- 
considered effect tha t could bias the meas ured age spread 
is that presented by Baraffe et al. (20091. Under certain 
conditions, if a star accretes a significant amount of ma- 
terial it cannot adjust its structure quickly enough, so 
appears underluminous for a star of its mass. This effect 
will have a memory, which may persist for several Myr 
according to the authors, and will cause some young stars 
to appear older than the are. This would also cause an 
additional apparent spread in the H-R diagram even in 
case of a coeval pop ulation. The scenario proposed by 
Baraffe et al. (2009), however, requires episodes of vig- 
orous mass accretion, with a high M. With no direct 
measurements of accretion in LH 95, and no information 
about the past history of accretion, as a function of mass, 
we are unable to estimate if this scenario could play any 
role in our region, and if so, to what extent. 

8. SUMMARY 

In this study we utilize the photometric catalog of 
sources in the young LMC stellar association LH 95, 
based on HST/ACS photometry in F555W and F814W 
filters (corresponding to standard V- and /-bands). The 
catalog, which covers stellar masses as low as 0.2 M Q was 
used in Paper I to derive for the first time the sub-solar 
IMF in the LMC. Here we focus on the age of the sys- 
tem and the determination of the age-spread from the 
CMD-broadening of its PMS stars. We first refine the 
decontamination of our sample from field stars by elim- 
inating any remaining field MS stars still present in the 
first catalog of true stellar members (Paper I) . Our pho- 
tometric catalog of member-stars of LH 95 covers ~ 2 000 
members in the pre-main-sequence (PMS) and upper- 
main-sequence (UMS) regions of the CMD. 

We develop a maxi mum-likelihood method, similar to 
the that described in (Naylor & Jeffries 2006) to derive 
the age of the system accounting simultaneously for pho- 
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tometric errors, unresolved binarity, differential extinc- 
tion, accretion, stellar variability and crowding. The ap- 
plication of the method is performed in several steps: 
First, a set of PMS isochrones is converted to 2D prob- 
ability distributions in the magnitude-magnitude plane; 
one model distribution is constructed for every consid- 
ered age. A well populated (N > 10 6 stars) synthetic 
population is thus created by applying the aforemen- 
tioned sources of displacement in the CMD - except for 
photometric errors - to the isochrone for each chosen 
age. A maximum-likelihood method is then applied to 
derive the probability, for each observed star to have a 
certain age, taking its photometric uncertainty into ac- 
count. The multiplication of the likelihood functions of 
all the stars provides an age probability function for the 
entire cluster, from where the most probable age of the 
cluster is derived. 

We apply our method to the stellar catalog of LH 95, 
assuming the IMF and reddening distribution for the sys- 
tem, as they were derived in Paper I. We assume distri- 
butions of the variations in both 771555 and 777,814 mag- 
nitudes caused by accretion and variability due to dark 
spots by utilizing data for the Orion Nebula Cluster from 
the literature. We model age-dependant probability dis- 
tributions using tw o sets of evolutionary models, that by 
Siess et al. (2000), for the metallicity of the LMC, and 
a new gri d computed with the FRANEC stellar evolu- 
tion c ode ( Chieffi fc Straniero|19 89; Dcgl'Innocenti et al 



2008[ ) for the same metallicity. We evaluate the accu- 
racy of our age estimation and we find that the method 
is fairly accurate in the PMS regime, while the precision 
of the measurement of the age is lower at higher lumi- 
nosities and old ages. Our method cannot constrain the 
age of UMS stars due to the small dependence of CMD 
position with age. Our treatment showed that the age 
determination is very sensitive to the assumed evolution- 
ary model and binary fraction; the age of LH 95 is found 
to vary from ~2.8 Myr to ~4.4 Myr, depending on these 
factors. Assuming / = 0.5, the average binary fraction 
in our mass range (Lada 2006), the best-fit age of the 



system derived for this binary fraction is determined to 
be of the order of 4 Myr. Despite the sparse and clumpy 
spatial distribution of PMS stars in LH 95, we find no 
evidence of systematic variations of age with location. 

With our method we are also able to assess the true 
appearance of an age-spread among the PMS stars of 
LH 95, as an indication that star formation in the system 
occurred not spontaneously but continuously. Our treat- 
ment allowed us to disentangle any real age-spread from 
the apparent CMD-spread caused by the physical and 
observational biases which dislocate the PMS stars from 
their theoretical positions in the CMD (see § |5.1|). We 



find that LH 95 hosts an overall age-spread well modeled 
by a Gaussian distribution with FWHM = 2.8-4 Myr. 

We measure separately ages and age spread as a func- 
tion of luminosity. For the intermediate masses (M ~ 
IMq) we find younger average ages and larger age spread, 
whereas for low- mass stars (M ~ 1M Q ), older ages and 
more modest age spreads. We assign this dependence of 
age with mass not to a real effect but to imperfections 
of the PMS evolutionary models, which tend to predict 
lower ages for the intermediate masses, and higher ages 
for low-mass stars. In any case, under the validity of 
the assumptions used in our modeling, we rule out the 
possibility that the population of LH 95 is coeval. 

Finally, we remark that, despite the method we 
have presented provides an exceptional tool for an 
unbiased investigation of ages and age spreads in young 
PMS regions, the results are highly dependent on the 
modeling of the assumed sources of apparent spread in 
the CMD. Only a precise characterization of each of 
these effects, in different regions, will enable us to give 
the final answer on the long-standing questions about 
the duration of star formation. 
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